Sin duda alguna los automóviles han sido y son parte esencial de la sociedad, su movilidad y su transporte, desafortunadamente los siniestros causados por ellos son parte inherente de su uso y el impacto de ellos repercute de manera directa o indirecta en la vida de los involucrados, pero también de las personas que los rodean.
De acuerdo con un reporte de la Organización Mundial de la Salud en 2017, los accidentes viales son causantes de más de 1.3 de millones de muertes al año, y un estimado de alrededor de 50 millones de personas con traumatismos no mortales siendo así la quinta causa de los fallecimientos mundiales anuales.
Aunado a lo anterior, la OMS ahonda en su reporte, indicando que el 93% de los accidentes de tránsito suceden en los países de ingresos bajos y medianos. De este 93%, 20% acaece en la región de las Américas, en donde las estadísticas reportan más de 155’000 muertes por año por accidentes automovilísticos.
Por su parte México y de acuerdo con los datos del Instituto Nacional de Salud Pública (INSP), el país se encuentra en el séptimo lugar a nivel mundial en el rubro de siniestros automovilísticos mortales, y ocupa el tercer lugar en la región de Latinoamérica en este mismo rubro
Además de las repercusiones sociales y humanas, los accidentes viales, representan en números un estimado de más de 220 mil millones de pesos en costos, cifra la cual supera el 2% del PIB del país. Eventualmente los costos no son todos absorbidos por las instituciones del gobierno, y es el propio usuario en quién recae la mayor cantidad de responsabilidad económica asociada al siniestro. De nuevo, este último factor es un atenuante para la sociedad puesto que, continuando con los estudios presentados por el INSP, en promedio el gasto por accidente alcanza un costo de alrededor de 25´000 pesos, y si a esto se le suma que derivado de los siniestros automovilísticos se suman alrededor de 50´000 personas al grupo de personas con discapacidad motora, de los cuáles 70% no vuelve a conseguir empleo. Los costos en seguridad social, reparación de daños, entre otros regresan para ser solventados por las instituciones de gobierno.
Bajo este tenor, y con el fin de atenuar la situación descrita anteriormente, resulta útil identificar cuáles son los factores que tienen mayor influencia en los choques automovilísticos, cuáles son los perfiles más propensos a generar un choque, y simultáneamente diseñar estrategias que permitan atacar las causas raíz de los siniestros automovilísticos en México.
import pandas as pd
import matplotlib.pyplot as plt
import seaborn as sns
import numpy as np
import os
import urllib.request
url="https://www.inegi.org.mx/contenidos/programas/accidentes/datosabiertos/atus_anual_csv.zip"
file_name="atus_anual_1997_2019"
def descarga(url,file_name,n_file=0):
import requests, zipfile, io,os.path
if not os.path.exists(file_name):
r = requests.get(url)
z = zipfile.ZipFile(io.BytesIO(r.content))
if n_file==0:
z.extractall("")
else:
z.extractall(file_name)
return print("Descarga de archivos completa")
descarga(url,file_name)
Lectura de los datos pertenecientes a los años 2010 - 2019
from os import listdir
filepaths = [f for f in listdir("atus_anual_1997_2019/conjunto_de_datos") if f.endswith('.csv')]
lista_archivos=list(map(lambda x: "atus_anual_1997_2019/conjunto_de_datos/"+x, sorted(filepaths,reverse=True)[:10]))
df = pd.concat(map(lambda archivo: pd.read_csv(archivo,index_col=False),lista_archivos))
Dataset Original perteneciente a accidentes en todo Mexico del 2010 al 2019
df.shape
Tras la conjunción de todos los archivos obtuvimos un dataframe con 45 campos y 3,841,493 registros
df.head()
Para la acotación correspondiente a nuestro proyecto, Filtramos los datos pertenecientes a la Ciudad de Mexico cuyo codigo de entidad es 9
df_cdmx = df[df.ID_ENTIDAD==9]
df_cdmx = df_cdmx.reset_index(drop=True)
df_cdmx.shape
Tras el filtro para seleccionar únicamente a la CDMX, se obtuvo un datframe más pequeño con 136,812 registros. Continuando con el análisis exploratorio de los datos, examinamos el nombre de las columnas de la siguiente forma
df_cdmx.columns
A continuación se muestra la descrpición de los campos, la cual nos ayuda a entender mejor el significado de los datos de cada columnna.
diccionario_datos=pd.read_csv("atus_anual_1997_2019/diccionario_de_datos/diccionario_de_datos_atus_anual_1997_2019.csv",
index_col="COLUMNA",usecols=[0, 1])
diccionario_datos.head()
Lista de los campos y su descripción.
def imprimir_descripcion(df, width):
with pd.option_context('display.max_colwidth', width):
print(df)
imprimir_descripcion(diccionario_datos, 800)
Igualmente de la siguiente forma podemos elegir el nombre de la columna de la cual queremos obtener la descripcion
list(diccionario_datos.loc["ID_MUNICIPIO"])
Examinando la Columna ID_MUNICIPIO podemos observar que son solo numeros, por lo cuel seria mejor trabajar con nombre de los municipios en lugar de numero. Por lo cual cambiaremos esta columna. Para lograrlo importamos la informacion perteneciente al nombre del municipio asociado con su ID.
#Lectura del nombre pertenecientes a las delegaciones
delegacion_nombre=pd.read_csv("atus_anual_1997_2019/catalogos/tc_municipio.csv",index_col=False)
Examinando el archivo podemos observar que son solo 3 columnas por lo cual filtramos los datos para obtener solo los municipios pertencecientes a la CDMX los cuales deben de ser 16
delegacion_nombre=delegacion_nombre[(delegacion_nombre.ID_ENTIDAD==9)&(delegacion_nombre.ID_MUNICIPIO!=999)]
delegacion_nombre.head()
Agregamos el nombre de las delegaciones al Dataframe que estamos manejando
df_cdmx=df_cdmx.merge(delegacion_nombre[['ID_MUNICIPIO', 'NOM_MUNICIPIO']], how='left', on='ID_MUNICIPIO')
df_cdmx.head()
Examinando las columnas observamos que existen columnas separadas para los elementos de pertenecientes a fecha del accidente. Es decir, existe una columna para año, otra para mes, dia, etc. Por lo cual seria adecuado combinar estas columnas para obtener la fecha completa en la cual el accidente tuvo lugar
df_cdmx = df_cdmx.reset_index(drop=True)
concatenar = lambda anio,mes,dia,hora,minu: str(anio)+"/"+str(mes)+"/"+str(dia)+" "+ str(hora)+":"+str(minu)
df_cdmx["FECHA"] = list(map(concatenar,df_cdmx.ANIO,df_cdmx.MES,df_cdmx.ID_DIA,df_cdmx.ID_HORA,df_cdmx.ID_MINUTO))
df_cdmx["FECHA"] = pd.to_datetime(df_cdmx["FECHA"])
df_cdmx.head(5)
Despues de realizar los pasos anteriores hemos agregado mas columnas a nuestro dataset, por lo cual empezaremos a revisarlas, para ver que no haya errores ni NaNs.
df_cdmx.info()
Observando el dataset podemos notar que no tenemos valores nulos o NaNs. De igual forma muchas columnas estan en formato object, por lo cual en su mayoria son texto. Analizando la documentacion perteneciente al Dataset, sabemos que la mayoria de columnas del tipo object pertenecen a datos del tipo categorico, por lo cual las opciones son limitadas. Examinando estas columnas podriamos encontrar algun error en el momento de su captura.
columnas_object=df_cdmx[['COBERTURA', 'DIASEMANA','URBANA','SUBURBANA','TIPACCID',
'CAUSAACCI','CAPAROD','SEXO','ALIENTO','CINTURON','CLASACC',
'ESTATUS','NOM_MUNICIPIO']]
Mediante la siguiente funcion podemos examinar con mas detalle los valores unicos de cada columna tipo object
def rstr(df): return df.apply(lambda x: [x.unique()])
imprimir_descripcion(rstr(columnas_object), 800)
Podemos observar que en la columna de dias, existen errores de captura ya que existe Sabado y Sábado, por lo cual corregiremos esto.
map_value={"Sábado":"Sabado",
"Miércoles":"Miercoles",
"lunes":"Lunes"
}
df_cdmx['DIASEMANA']=df_cdmx['DIASEMANA'].map(map_value).fillna(df_cdmx['DIASEMANA'])
df_cdmx['DIASEMANA'].unique()
Una vez analizadas las columnas object, analizaremos las columnas numericas mediante la instruccion describe
df_cdmx.describe()
Analizando la informacion podemos observar que las columnas pertenecientes a numero de muertos (NEMUERTO) y numero de heridos (NEHERIDO) estan llenas de ceros, por lo cual se podrian depurar del set de datos.
De igual forma, las columnas PASAHERIDO,PEATMUERTO,PEATHERIDO,CICLMUERTO,CICLHERIDO,OTROMUERTO,OTROHERIDO estan en su mayoria llenas de ceros de acuerdo a lo que observamos ya que el 75% de los datos (tercer quartil) son 0. Podria ser buena idea depurarlos, sin embargo podrian ser de utilidad para analizis posteriores, por lo cual se conservaran
df_cdmx=df_cdmx.drop(['NEMUERTO',"NEHERIDO"], axis=1)
#generacion del archivo de salida
#compression_opts = dict(method='zip',archive_name='Accidentes_vehiculares_cdmx.csv')
#df_cdmx.to_csv("Accidentes_vehiculares_cdmx.zip",compression=compression_opts)
Una vez arreglados los datos podemos empezar a realizar algunas preguntas referentes a ellos. Para este proceso se realizo un Brainstorming con el fin de registrar la mayor cantidad de preguntas derivadas de la observación de nuestros datos, este ejercicio permitió empatar preguntas claves por grupos y tras la evualación final, se llego a la obtención de las preguntas que veremos a continuación.
¿Que delegacion/municipio resgistro mas accidentes?
Podemos conocer esto haciendo una agrupacion por Delegacion y en base al año y ordenandola de manera descendiente en base a la delegacion y el año
delagacion_accidentes=df_cdmx.groupby("NOM_MUNICIPIO")["ANIO"].value_counts()
delagacion_accidentes.sort_index(level=['NOM_MUNICIPIO','ANIO'],inplace=True)
delagacion_accidentes.head(30)
Resulta un poco dificil de analizar los datos, podemos graficarlos para analizarlos mejor de la siguiente manera
# Librearia necesaria para usar ggplor en python
#!pip install plotnine
from plotnine import *
%matplotlib inline
df_gra=df_cdmx.groupby(["NOM_MUNICIPIO","ANIO"]).size()#["FECHA_DIA"].value_counts()
df_gra=df_gra.reset_index()
df_gra.columns =["Municipio","Año","count"]
df_gra
plt.style.use('seaborn')
sns.set(font_scale=1.5)
fig = plt.figure(figsize=(12, 8))
ax = fig.add_subplot()
sns.barplot(x="Municipio",y="count", ax=ax, palette='husl',hue="Año",data=df_gra);
ax.set_ylabel('Numero de Accindentes')
ax.set_title('Histograma de Accidentes por Delegacion y año')
ax.set_xticklabels(ax.get_xticklabels(), rotation=85)
plt.legend(title="Año",bbox_to_anchor=(1.05, 1), loc='upper left', borderaxespad=0.);
De manera grafica se puede analizar un poco mejor los daTos y podemos observar las delegaciones con mas accidentes, las cuales son las delegaciones Cuauhtémoc,Iztapalapa y Gustavo A. Madero. Para tener un panorama mas claro podemos hacer otros analis de los datos, por ejemplo realizar un historgrama de los accidentes totales del 2010 al 2019 de acuerdo a la delegacion.
df_out=df_cdmx.groupby("NOM_MUNICIPIO").size().sort_values(ascending=False)
#df_out.to_csv("df.csv")
df_graph=df_cdmx.groupby("NOM_MUNICIPIO").size().sort_values(ascending=False)
sns.set_style('whitegrid')
fig = plt.figure(figsize=(8, 8))
sns.set(font_scale=1.2)
sns.set_style('white')
ax = sns.barplot(df_graph.values,df_graph.index)
ax.set_title("Histograma de Accidentes totales por Delegacion")
ax.set(xlabel="Accidentes")
ax.set(ylabel="Municipio");
Adicionalmente podemos realizar un mapa cloroplético para visualizar mejor el porcentaje de accidentes totales entre el año 2010 al 2019
!pip install folium
import folium
df_graph=df_graph.reset_index()
df_graph.columns =["Municipio","count"]
Para realizar el mapa cloropletico es necesario conseguir los datos GeoJson, los cuales son datos para el mapeo de las delegaciones, es decir son datos geograficos. Estos datos se encuentran en el dataset MX-CMX.json el cual esta en la misma carpeta que este Jupyter Notebook
#El archivo MX-CMX.json contiene los poligonos necesarios para graficar las delegaciones
state_data=df_graph
state_data.loc[7,'Municipio']="Alvaro Obregón"
state_data.loc[14,'Municipio']="Magdalena Contreras"
state_data
import json
#Revisando el GeoJson data para encontrar en que llave esta guardado el nombre de la delegacion
with open("MX-CMX.json",encoding="utf8") as f:
#map_data = f.readlines()
s = f.read()
s = s.replace('\t','')
s = s.replace('\n','')
s = s.replace(',}','}')
s = s.replace(',]',']')
data = json.loads(s)
data['features'][0]["properties"]["NAME_2"]
#Creando una lista con el nombre de los municipios extraidos en los GeoJson en el orden que aparecen
json_map_file = []
for i in range(len(data['features'])):
json_map_file.append(data['features'][i]['properties']['NAME_2'])
json_map_file = pd.DataFrame({'Municipio': json_map_file})
#Combinando los datos, para poder obtener los valores del numero de accidentes en el orden del GeoJson
mun_mapping = json_map_file.merge(state_data, on='Municipio')
#mun_mapping
#Creando una columna de datos que contenga el nombre del municipio asi como el porcentaje de accidentes.
mun_mapping["Percent"]=mun_mapping["Municipio"]+": "+((mun_mapping['count']/mun_mapping['count'].sum()) * 100).map(lambda x:str(round(x,2))+"%")
#[data['features'][idx]['properties']['tooltip1']=tooltip_text[idx] for x in mun_mapping["Percent"]]
for idx in range(len(mun_mapping["Percent"])):
data['features'][idx]['properties']['tooltip1'] = mun_mapping["Percent"][idx]
# map_data['features'][idx]['properties']['tooltip2'] = tooltip_text2[idx]
#Verificando que se haya creando el valor de la llave en el archivo Json
print(data['features'][0]['properties'])
#Grafico Iteractivo
bins = list(state_data["count"].quantile([0, 0.25, 0.5, 0.75, 1]))
m = folium.Map(location=[19.32847, -99.12766],tiles='cartodbpositron', zoom_start=10,width='80%', height='80%')
folium.Choropleth(
geo_data=data,
name='choropleth',
data=state_data,
columns=['Municipio', 'count'],
key_on='feature.id',
fill_color='BuPu',
fill_opacity=0.7,
line_opacity=0.5,
legend_name='Numero de accidentes',
bins=bins,
reset=True,
).add_to(m)
#folium.Marker(
# location=[19.5040652077,-99.1158642087],
# popup="Mt. Hood Meadows",
# icon=folium.Icon(icon="cloud"),
#).add_to(m)
folium.LayerControl().add_to(m)
style_function = lambda x: {'fillColor': '#ffffff',
'color':'#000000',
'fillOpacity': 0.1,
'weight': 0.1}
highlight_function = lambda x: {'fillColor': '#000000',
'color':'#000000',
'fillOpacity': 0.50,
'weight': 0.1}
NIL = folium.features.GeoJson(
data,
style_function=style_function,
control=False,
highlight_function=highlight_function,
tooltip=folium.features.GeoJsonTooltip(['tooltip1'] ,labels=False)
)
m.add_child(NIL)
m.keep_in_front(NIL)
folium.LayerControl().add_to(m)
m
De acuerdo a lo obtenido podemos comprobar las delegaciones con mas accidentes, en especial la delegacion Cuauhtemoc la cual registra los mayores accidentes. Seguido de ella se encuentra la delegación de Iztapalapa y Gustavo A. Madero. Por su parte las delegaciones de Cuajimalpa de Morelos, La Magdalena Contreras y Milpa Alta son las delegaciones con menos choques.
Tras realizar esta gráfica, concluimos que sería interesante también comparar la cantidad de vehiculos y población por delegación, con el fin de ahondar en el análisis.
¿Cual es el promedio de accidentes por delegacion del año 2010 a 2019?
Seria Interesante saber los accidentes promedio de acuerdo a cada delegacion, por lo cual realizaremos lo siguiente mediante la funcion groupby
av_e=df_cdmx.groupby("NOM_MUNICIPIO")["ANIO"].value_counts()
av_e.groupby("NOM_MUNICIPIO").mean().sort_values(ascending=False)
Como podemos ver el Promedio de accidentes de las delegaciones Cuauhtémoc,Iztapalapa y Gustavo A. Madero son las que tienen mayores promedios
¿Cuales son las tres delegaciones con mas accidentes y como es su tendencia?
Como lo hemos visto en las secciones previas, las delegaciones Cuauhtémoc,Iztapalapa y Gustavo A. Madero son las que registran mas accidentes. Ahora nos gustaria saber como es su tendencia a lo largo del periodo 2010 al 2019. Primero empezaremos ultizando el dataset creado en las secciones anteriores
delagacion_accidentes=delagacion_accidentes.loc[["Cuauhtémoc","Iztapalapa","Gustavo A. Madero"]]
delagacion_accidentes=delagacion_accidentes.reset_index(name="Num_Acc")
delagacion_accidentes
De manera grafica podriamos observar los resultados obtenidos y con el fin de facilitar la visualización de las tendencias se realizo un gráfico referente a ellos.
plt.style.use('ggplot')
fig = plt.figure(figsize=(10, 10))
sns.set(font_scale=1.2)
ax = sns.barplot(data=delagacion_accidentes, x="ANIO", y="Num_Acc", hue="NOM_MUNICIPIO")
ax.set_title('Delegaciones con mas accidentes', fontsize=20, pad=15)
ax.set(xlabel='Año', ylabel='Numero de Accidentes')
new_title = 'My title'
plt.legend(title='Delgacion');
(ggplot(delagacion_accidentes, aes(x="factor(ANIO)",y="Num_Acc", group = 1))
+geom_point()+geom_smooth(method='lm',span=.2, se=False)+facet_wrap("NOM_MUNICIPIO",ncol=1)
+ylab("Numero de Accidentes") + xlab("Año")
)
Como podemos observar las tendencias para las tres delegaciones es a la baja, por lo cual podriamos suponer que los accidentes en estas 3 delegaciones disminuiran en los proximos años. Durante los años 2012, y 2013 podemos observar que alcanzaron valores cercanos a los 2000 accidentes y en el 2019 únicamente la delegación Cuauhtémoc supero los 1500 accidentes.
¿En que año se han registrado mas accidentes en la CDMX y cual es su tendencia?
años_accidentes=df_cdmx.groupby("ANIO").size()
años_accidentes.sort_values(ascending=False)
Como podemos observar el año en el cual se registraron mas accidentes fue el 2012, sin embargo seria interesante ver las tendencias que tienen los accidentes, para esto haremos la siguiente grafica
años_accidentes=años_accidentes.reset_index(name="Num_Acc")
años_accidentes
años_accidentes.corr(method="pearson")
Para este caso específico se decidió analizar con el coeficiente de Pearson la relación entres los años y el número de accidentes registrados. Este podría ser un análisis preeliminar sobre la tendencia en el número de choques por año; de acuerdo al valor obtenido -0.87 podemos decir de primera instancia que a medida que aumenten los años, la cantidad de accidentes irán a la baja, sin embargo sería muy arriesgado sacar conclusiones a partir de este único indicador, por lo que en análisis posteriores, se complementará con el cooeficiente de spearman, la ecuación del gráfico y/o análisis más robustos como series de tiempos.
#plt.style.use('seaborn')
#fig = plt.figure(figsize=(5, 5))
#ax = fig.add_subplot()
#fig, ax = plt.subplots()
#ax2 = ax.twinx()
#ax.barplot(data=años_accidentes, x="ANIO", y="Num_Acc")
#ax2.regplot(x='ANIO', y='Num_Acc', data=años_accidentes, ax=axs[0])
#plt1=ax.barplot(data=años_accidentes, x="ANIO", y="Num_Acc");
#plt1 = ax.bar(value_counts.index, value_counts['M'], label='M',
# color=["#7788AA","#4E638E","#2E4372","#152A55"])
#plt2 = ax.bar(value_counts.index, value_counts['F'], bottom=value_counts['M'],
# color=["#FFD0AA", "#D4996A", "#AA6B39", "#804415"])
#ax.set_ylabel('count')
#ax.set_title('Conteo de frecuencia de 4 deportes Olímpicos', fontsize=13, pad=15);
#plt.legend((plt1[0], plt2[0]), ('Men', 'Women'));
#ax.set_ylim(0, 4500);
#sns.barplot(data=años_accidentes, x="ANIO", y="Num_Acc");
#sns.regplot(x='ANIO', y='Num_Acc', data=años_accidentes, ax=ax)
#sns.lineplot(x="ANIO", y="Num_Acc", data=años_accidentes);
(
ggplot(años_accidentes) +
aes(x="factor(ANIO)",y="Num_Acc", group = 1) + geom_col(alpha= 0.2, fill = "blue")+
ggtitle("Histograma de Accidentes en la CDMX") +geom_point()+geom_smooth(method='lm',span=.2)+
ylab("Accidentes totales") +
xlab("Años")
)
Similar a lo observado en la tendecia de las gráficas de las delegaciones con mayor cantidad de accidentes registrados y en concordancia con la relación del coeficiente de Pearson, se observa en los datos que la tendencia de los accidentes en la CDMX tiene una tendencia a disminuir, esto se debera validar posteriormente con modelos de predicción en futuros análisis
¿Quiénes son mas propensos a sufrir accidentes, los hombres o las mujeres y como es la tendencia?
Seria interesante analizar quienes son mas propensos a sufrir accidentes, para esto empezaremos analizando la cantidad de accidentes dependiendo del sexo.
df_cdmx.groupby('SEXO').size()#.map(lambda x:x*100/df_cdmx["ID_EDAD"].count()).loc[[0,99]]
A simple vista se puede observar que la mayoria de accidentes son causados por hombre. De igual forma, existe una gran cantidad de accidentes en los cuales el conducto se da a la fuga, Incluso superando los registros de mujeres involucradas en accidentes. Podriamos calcular el porcentaje de cada una de estas categorias para poder analizar mejor los datos.
df_cdmx.groupby('SEXO').size().apply(lambda x:x*100/df_cdmx["SEXO"].count())
Tal como lo habias observado el porcentaje de Mujeres involucradas en un accidente es mucho menor al de hombres, al igual que el de conductores que se dieron a la fuga. El siguiente paso es analizar las tendencias de accidentes tanto para hombres como mujeres, para saber si van a la baja o la alta.
filtro_sexo_accidentes=df_cdmx.groupby("SEXO")["ANIO"].value_counts().loc[["Hombre","Mujer"]]
filtro_sexo_accidentes=filtro_sexo_accidentes.reset_index(name="Num_Acc")
filtro_sexo_accidentes
(ggplot(filtro_sexo_accidentes, aes(x="factor(ANIO)",y="Num_Acc",group = 1))
+geom_col(alpha = 0.3)+facet_grid('SEXO ~ .',scales="free")+geom_point()+geom_smooth(span=.2,se=False)
+ylab("Numero de Accidentes") + xlab("Año")
)
Adicionalmente podemos analizar el caso de las personas que se fugaron y su tendencia.
filtro_fuga_accidentes=df_cdmx.groupby("SEXO")["ANIO"].value_counts().loc["Se fugó"]
filtro_fuga_accidentes=filtro_fuga_accidentes.reset_index(name="Num_Acc")
(
ggplot(filtro_fuga_accidentes) +
aes(x="factor(ANIO)",y="Num_Acc", group = 1) + geom_col(alpha= 0.2, fill = "blue")+
ggtitle("Tedencia de la fuga de personas en accidentes en la CDMX") +geom_point()+geom_smooth(method='lm',span=.2)+
ylab("Accidentes totales") +
xlab("Años")
)
De acuerdo a los graficos podemos observar una tendencia clara de disminucion de accidentes en el caso de hombres. En el caso de las mujeres la tendencia no es tan clara ya que hubo un aumento considerable de accidentes en el 2013, y en años posteriores la tendencia fue ligeramente a la baja. Por lo anterior podriamos suponer unas tendencias a la baja de los dos sexos en cuanto a los accidentes.
En cuanto a la tendencia de las personas a fugarse de los accidentes, los datos indican en general una tendencia de dismonucion. Sin embargo en los ultimos años el numero de personas que se fugan de los accidentes ha comnezado a aumentar ligeramente. Pero en general podriamos suponer que va a la baja.
¿Cuáles son los meses con mayores accidentes a los largo de los años?
Se realizó un análisis respecto a los meses con mayor número de accidentes a lo largo de los 10 años observados.
Para empezar de agruparon los datos del dataframe original por año y mes, posteriomente se concatenaron en un nuevo dataframe con únicamente los máximos de accidentes de todos los años, y se modificaron los meses para una mayor comprensión de los datos.
df_anios = df_cdmx.groupby(['ANIO','MES']).size()
print(df_anios)
df_maximos_anios = pd.DataFrame()
for anio in range(df_cdmx['ANIO'].unique().max(),df_cdmx['ANIO'].unique().min()-1,-1):
df_maximos_anios = pd.concat([df_maximos_anios,df_anios[anio][df_anios[anio] == df_anios[anio].max()]])
df_maximos_anios['MES'] = df_maximos_anios.index
int_a_mes = {1 :'Enero',2 :'Febrero',3 : 'Marzo',4 : 'Abril', 5 : 'Mayo', 6 : 'Junio',
7 : 'Julio', 8 : 'Agosto', 9 : 'Septiembre',10 : 'Octubre', 11 : 'Noviembre', 12 : 'Diciembre'}
df_maximos_anios['MES'] = df_maximos_anios['MES'].map(int_a_mes)
df_maximos_anios.set_index(df_cdmx['ANIO'].unique(),inplace=True)
df_maximos_anios.rename(columns={0 : 'Max_acci'},inplace=True)
df_maximos_anios.index.name = 'ANIO'
df_maximos_anios
En el código aqui presente se obtuvieron los meses con mayor número de accidentes en los 10 años y se observó que el mes con más accidentes fue Obtubre en 6 años del total observado, seguido por Agosto, Septiembre y Junio respectivamente.
df_maximos_anios.groupby('MES').size().sort_values(ascending=False)
df_gra=df_cdmx.groupby(["ANIO","MES"]).size()#["FECHA_DIA"].value_counts()
df_gra=df_gra.reset_index()
df_gra.columns =["Año","Mes","count"]
plt.style.use('seaborn')
sns.set(font_scale=1.5)
fig = plt.figure(figsize=(12, 8))
ax = fig.add_subplot()
sns.barplot(x="Mes",y="count", ax=ax, palette='Paired',hue="Año",data=df_gra);
ax.set_ylabel('Numero de Accindentes')
ax.set_title('Histograma de Accidentes por mes y año')
labels=["Ene","Feb","Mar","Abr","May","Jun","Jul","Ago","Sep","Oct","Nov","Dic"]
ax.set_xticklabels(labels)
plt.legend(title="Año",bbox_to_anchor=(1.05, 1), loc='upper left', borderaxespad=0.);
#ax = sns.barplot(x="day", y="total_bill", hue="sex", data=tips)
Se grafican los datos de los accidentes contra los meses segun los años observados para comprobar los resultados arrojados en el proceso anterior.
¿Cuáles son los dias de la semana con mayores accidentes a los largo de los años?
Se realizó un analisis respecto a los dias con mayor número de accidentes a lo largo de los 10 años observados.
df_cdmx["FECHA_DIA"]=df_cdmx["FECHA"].dt.weekday
#df_cdmx.groupby("FECHA_DIA").size()
df_gra=df_cdmx.groupby(["FECHA_DIA","ID_HORA"]).size()#["FECHA_DIA"].value_counts()
df_gra=df_gra.reset_index()
df_gra.columns =["Dia","Hora","count"]
df_gra
plt.style.use('ggplot')
g = sns.FacetGrid(data = df_gra, col = 'Dia', col_wrap = 3)
g.map(plt.bar, 'Hora','count')#.set(yscale = 'log')
g.set_xlabels('Hora del Dia')
g.set_ylabels('Numero de Accidentes')
days=["Lunes","Martes","miercoles","Jueves","Viernes","sabado","domingo"]
#ax.set(xlabel='Año', ylabel='Numero de Accidentes');
for ax in g.axes.flatten():
ax.tick_params(labelbottom=True)
#ax.set(xlabel='Hora del Dia')
axes = g.axes.flatten()
[axes[x].set_title(days[x]) for x in range(7)];
En esta gráfica observamos la cantidad de accidentes a lo largo del día y evaluado por día de la semana. El Lunes es nuestro día 0 y el Domingo el día 6. Podemos observar que la cantidad de accidentes se comporta de manera ascendente conforme transcurre el día con excepción de horas del día específicos como 9 de la mañana. El Viernes podemos ver valores más altos y mantenidos, sin embargo Míercoles y Jueves destacan por tener los valores más altos registrados.
Con los resultados de las gráficas anteriores se procedió a agrupar los datos del dataframe original por año y dia, posteriomente se concatenaron en un nuevo dataframe con únicamente los máximos de accidentes de todos los años.
df_anios_dias = df_cdmx.groupby(['ANIO','DIASEMANA']).size()
df_maximos_anios_dias = pd.DataFrame()
for anio in range(df_cdmx['ANIO'].unique().max(),df_cdmx['ANIO'].unique().min()-1,-1):
df_maximos_anios_dias = pd.concat([df_maximos_anios_dias,df_anios_dias[anio][df_anios_dias[anio] == df_anios_dias[anio].max()]])
df_maximos_anios_dias
A primera vista, se observa que a lo largo de los años observados el dia viernes destaca como el dia donde hay mayores accidentes en cada una de las observaciones.
df_maximos_anios_dias['DIAS'] = df_maximos_anios_dias.index
df_maximos_anios_dias.set_index(df_cdmx['ANIO'].unique(),inplace=True)
df_maximos_anios_dias.rename(columns={0 : 'Max_acci'},inplace=True)
df_maximos_anios_dias.index.name = 'ANIO'
df_maximos_anios_dias
Una vez con los datos acomodados correctamente se procedió a obtener la cantidad de años en donde cada día tuvo el máximo de accidentes, siendo el dia viernes con mayor número de accidentes ya que en 9 años del total observado este dia destaca con mayores accidentes, seguido por el dia jueves.
df_maximos_anios_dias.groupby('DIAS').size().sort_values(ascending=False)
¿Cual es la edad mas propensa a sufrir accidentes?
Para responder esta pregunta, primero tenemos que contabilizar el numero de accidentes en base a la edad, del año 2010 al 2019, depues ordenaremos los datos de manera descendiente para poder analizar las edades con mas accidentes
df_cdmx.groupby("ID_EDAD")["ID_EDAD"].count().sort_values(ascending=False)
De manera grafica podemos verlo de la siguiente manera
fig = plt.figure(figsize=(12, 8))
sns.set(font_scale=1.5)
ax=sns.distplot(df_cdmx['ID_EDAD'], kde=False, norm_hist=False, bins=100)
ax.set_title('Histograma de Accidentes en la CDMX', fontsize=20, pad=15)
ax.set(xlabel='Edad', ylabel='Numero de Accidentes');
Como podemos observar existe una gran cantidad de datos pertenecientes a la edad de 0 y 99 lo cual no hace sentido, esta gráfica es como se observan los datos sin ninguna manipulación, sin embargo se asume que son datos que no lograron recabarse de manera correcta por lo cuál se decidio abordar la situación de la manera mostrada abajo.
¿Que porcentaje de datos de edades de 99 y 100 años hay en los datos?
Apartir de la utima grafica obtenida de numero de accidentes en funcion de la edad se comienza a analizar los datos. Primero se observa que los valores de edad 0 y 99 albergan demasiados datos, de acuerdo a los registros del INEGI cuando el conductor se fugaba, los registros asignan el valor 0 a la edad. En cuanto a 99 significa que se desconoce la edad del conductor.
df_cdmx.groupby('ID_EDAD').size().map(lambda x:x*100/df_cdmx["ID_EDAD"].count()).loc[[0,99]]
Como podemos notar el porcentaje de datos pertenecientes a estas dos categorias es arriba del 30% de los datos por lo cual no es muy recomendable eliminar todas las filas pertenecientes a estos datos para futuras predicciones o algoritmos de prediccion. Sin embargo, podriamos analizar un poco los datos de edades excluyedo estas dos categorias
filtro_no_cero=df_cdmx["ID_EDAD"]!=0
filtro_no_noventa=df_cdmx["ID_EDAD"]!=99
df_cdmx_filtrado=df_cdmx[filtro_no_cero & filtro_no_noventa]
df_cdmx_filtrado.head()
fig = plt.figure(figsize=(12, 8))
sns.set(font_scale=1.5)
ax=sns.distplot(df_cdmx_filtrado['ID_EDAD'], kde=False, norm_hist=False, bins=90)
ax.set_title('Histograma de Accidentes en la CDMX', fontsize=20, pad=15)
ax.set(xlabel='Edad', ylabel='Numero de Accidentes');
Filtrando los datos, podemos observar un mejor comportamiento en la grafica, de igual forma podemos determinar la edad que ha sufrido mas accidentes (moda) y el promedio de edad
#df_cdmx_filtrado.groupby("ID_EDAD")["ID_EDAD"].count().loc[]
#Promedio de edad
df_cdmx_filtrado["ID_EDAD"].mean()
#Edad con mas Accidentes
df_cdmx_filtrado["ID_EDAD"].mode()
#top 5 edades
df_cdmx_filtrado.groupby("ID_EDAD")["ID_EDAD"].count().sort_values(ascending=False).head(5)
De acuerdo a los registros 30 años es la edad que tiene mas accidentes acumulados del 2010 al 2019, de igual forma podemos ver el top 5 de edades con mas accidentes acumulados siendo todas edades arriba de los 27 años.
De nuevo y en un posible anaálsis futuro sería recomendable evaluar la cantidad de accidentes por edad y la cantidad de población de esa edad que transita. Esto con el fin de eliminar sesgos, similar a lo que ocurre en el sexo.
¿Cual es la Tedencia de las edades mas propensas a sufrir accidentes del 2010 al 2019?
Para esto analizaremos las 3 edades con mas registros de accidentes y veremos como han variado a lo largo de los años
Tendencia_años=df_cdmx_filtrado.groupby('ID_EDAD')['ANIO'].value_counts()
Tendencia_30_años=Tendencia_años.loc[30].sort_index(ascending=True)
Tendencia_30_años
Tendencia_30_años.plot(xlabel='Año', ylabel='Numero de Accidentes');
T_30=Tendencia_30_años.to_frame()#.reset_index(level='ANIO',col_level=1, col_fill='genus')
T_30.columns = ['NUM_ACC']
T_30.reset_index(inplace=True)
T_30.corr(method="pearson")
Como se puede obsevar en la gráfica el comportamiento es un poco irregular. Definitavamente requerirá se lleve a cabo un análisis mucho más completo. Nuestro análisis preeliminar utilizando coeficiente de Pearson, tampoco nos permite ver alguna relación directa o indirecta significativa. Un ajuste polinomial de grado 5 probablemente podrá ajustarse a esta gráfica, aunque sería bueno poder realizar algun estudio con serie de tiempo para obtener resultados mas precisos.
Tendencia_40_años=Tendencia_años.loc[40].sort_index(ascending=True)
Tendencia_40_años
Tendencia_40_años.plot(xlabel='Año', ylabel='Numero de Accidentes');
Observamos en esta gráfica, que la tendencia de este grupo de edad va a la baja de acuerdo a los datos obtenidos en el 2019, al igual que los grupos de años anteriores en el 2013 fue cuando se obtuvo el valor máximo registrado de número de accidentes. Sin emabrgo, el coeficiente de Pearson es positivo pero al ser tan cercano a cero indica una correlación casi nula que es lógica, dado que pese a que son pocos datos hay mucha variación entre los años.
T_40=Tendencia_40_años.to_frame()#.reset_index(level='ANIO',col_level=1, col_fill='genus')
T_40.columns = ['NUM_ACC']
T_40.reset_index(inplace=True)
T_40.corr(method="pearson")
Tendencia_28_años=Tendencia_años.loc[28].sort_index(ascending=True)
Tendencia_28_años
Tendencia_28_años.plot(xlabel='Año', ylabel='Numero de Accidentes');
T_28=Tendencia_28_años.to_frame()#.reset_index(level='ANIO',col_level=1, col_fill='genus')
T_28.columns = ['NUM_ACC']
T_28.reset_index(inplace=True)
T_28.corr(method="pearson")
Observamos un comportamiento similar a las anteriores en los primeros años de la gráfica, con una tendencia a la baja en el último año, pero un 2017 en compración con los registros de los demás años, de nuevo es difícil genera una interpretación o predicción con el gráfico y el coeficiente de Pearson, definitavamente el análisis inferencial será de gran utilidad para poder dar una predicción más certera del comportamiento de este tipo de siniestro en los años venideros.
En base a los resultados obtenidos podemos obervar que para las edades de 40 y 28 años las tendencias parecen ir a la baja lo que nos podria indicar una reduccion en los accidentes a esta edad. En cambio para la edad con mas registros de accidentes que es 30 años, la tendencia parece indicar que seguira siendo la edad con mas accidentes.
¿Que vehiculo es el que causa mas accidentes?
Como se observo cuando se estaba analizando el Dataset, cada tipo de vehiculo tiene su propia columna por lo cual es necesario tranformar estas columnas en filas para poder realizar algunas agregaciones y calculos.
Para realizar esto, primero necesitamos elegir mas columnas que queremos trasnformar, las cuales se almacenaran en modelo_columnas.
modelo_columnas=["AUTOMOVIL","CAMPASAJ","MICROBUS","PASCAMION","OMNIBUS","TRANVIA",
"CAMIONETA","CAMION","TRACTOR","FERROCARRI","MOTOCICLET","BICICLETA","OTROVEHIC"]
columnas_filtradas=[i for i in list(df_cdmx.columns) if i not in modelo_columnas]
df_tipo_vehiculo=df_cdmx.melt(id_vars=columnas_filtradas,var_name ='TIPO_VEHICULO', value_name ='NUM_ACC')
df_tipo_vehiculo.head()
df_vehiculo_agg=df_tipo_vehiculo.groupby(["ANIO","TIPO_VEHICULO"])["NUM_ACC"].sum()
df_vehiculo_agg
Arreglando la serie para convertirla en un data frame, podemos graficarlo para poder visulaizar los resultados
df_vehiculo_agg=df_vehiculo_agg.to_frame().reset_index()
df_vehiculo_agg
plt.style.use('ggplot')
sns.set(font_scale=1.2)
fig = plt.figure(figsize=(12, 8))
ax = fig.add_subplot()
sns.barplot(x="NUM_ACC",y="TIPO_VEHICULO", ax=ax, palette='viridis',hue="ANIO",data=df_vehiculo_agg);
ax.set_ylabel('Tipo de vehiculo')
ax.set_xlabel('Numero de Accidentes')
ax.set_title('Histograma de Accidentes por tipo de vehiculo')
plt.legend(title="Año",bbox_to_anchor=(1.05, 1), loc='upper left', borderaxespad=0.);
#ax = sns.barplot(x="day", y="total_bill", hue="sex", data=tips)
!pip install plotly
Podemos visualizar los datos de una manera mejor utilizando Treemaps.
df_vehiculo_agg
df_vehiculo_agg["Tipo"]="Accidentes"
import plotly.express as px
fig = px.treemap(df_vehiculo_agg, path=["Tipo",'ANIO', 'TIPO_VEHICULO'], values='NUM_ACC')
fig.show()
En el grafico y analisis de datos podemos observar como, los accidentes de autos son los mas comunes. Lo cual es algo que ya esperabamos, ya que existen muchos mas automoviles en circulacion que los otros tipos de vehiculos.
df_vehiculo_agg.groupby("TIPO_VEHICULO")["NUM_ACC"].sum().sort_values(ascending=False)
Tambien calculando el acumulado de los accidentes, podemos determinar los vehiculos que mas accidentes generaron en este caso el top 3 de tios de vehiculos que causaron mas accidentes son Automovil, Camioneta de pasajeros y Motocicletas
¿Cual es la tendencia de los 6 tipos de vehiculos que causan mas accidentes?
Una vez que tenemos claro los tipos de vehiculos que causan mas accidentes, seria interesante observar si existe una tendencia en los tipos de vehiculos con mas accidentes y si esta es al aumento o disminucion
Primero empezaremos filtrando los 3 tipos de vehiculos con mas accidentes los cuales son Automovil(AUTOMOVIL), Camioneta de pasajeros (CAMPASAJ) y Motocicletas (MOTOCICLET)
filtro_top=df_vehiculo_agg.query('TIPO_VEHICULO=="AUTOMOVIL" or TIPO_VEHICULO=="CAMPASAJ" or TIPO_VEHICULO=="MOTOCICLET"')
(ggplot(filtro_top, aes("factor(ANIO)","NUM_ACC",group=1))
+geom_col(alpha= 0.2, fill = "blue")
+geom_point()+geom_smooth(method='lm',span=.2)
+facet_grid("TIPO_VEHICULO~ .",scales = 'free')
+theme(axis_text_x = element_text(angle=45, hjust=1))+ ggtitle("Accidentes por tipo de vehiculo y año") +
xlab("Año") +labs(fill = "YEAR")+ylab("Numero de accidentes")
)
A diferencia de la tendencia del camión de pasajeros y de los automóviles que presentan ambos una tendecia a la baja, las motocicletas tienen una ligera tendencia a la alta, especialmente en el año 2019.
De nuevo será útil valorar y analizar estas tendencias con la cantidad de parque vehicular relativo de cada tipo de vehículo.
Para analizar los 3 vehiculos pertenecistes al top 6 se realizar lo siguiente:
filtro_top=df_vehiculo_agg.query('TIPO_VEHICULO=="PASCAMION" or TIPO_VEHICULO=="CAMIONETA" or TIPO_VEHICULO=="MICROBUS"')
(ggplot(filtro_top, aes("factor(ANIO)","NUM_ACC",group=1))
+geom_col(alpha= 0.2, fill = "blue")
+geom_point()+geom_smooth(method='lm',span=.2)
+facet_grid("TIPO_VEHICULO~ .",scales = 'free')
+theme(axis_text_x = element_text(angle=45, hjust=1))+ ggtitle("Accidentes por tipo de vehiculo y año") +
xlab("Año") +labs(fill = "YEAR")+ylab("Numero de accidentes")
)
A diferencia del caso anterior, las tres gráficas de los tipos de automóviles presentan una línea de tendencia que va a la baja, para el caso del Microbus es espcialmente marcada esta tendencia, por lo cual podriamos de igual manera analizar que factores pudieron influir para que esto sucediera.
url= "https://www.inegi.org.mx/contenidos/programas/vehiculosmotor/datosabiertos/vmrc_anual_csv.zip"
file = "vmrc_anual_csv.zip"
if not os.path.isfile(file):
print("Descargando archivo")
r = urllib.request.urlopen(url)
f = open(file,"wb")
f.write(r.read())
f.close()
print("Se ha descargado el archivo")
else:
print("Ya existe el archivo")
if not os.path.exists("vmrc_anual_csv"):
os.mkdir("vmrc_anual_csv")
zf = zipfile.ZipFile(file, "r")
zf.extractall(path = "vmrc_anual_csv")
f=pd.read_csv("vmrc_anual_csv\conjunto_de_datos/vmrc_2019.csv",index_col=False)
filepaths = [f for f in listdir("vmrc_anual_csv\conjunto_de_datos") if f.endswith('.csv')]
lista_archivos=list(map(lambda x: "vmrc_anual_csv\conjunto_de_datos/"+x, sorted(filepaths,reverse=True)[:10]))
df = pd.concat(map(lambda archivo: pd.read_csv(archivo,index_col=False),lista_archivos))
Para la acotación correspondiente a nuestro proyecto, Filtramos los datos pertenecientes a la Ciudad de Mexico cuyo código de entidad es 9
df_cdmx_parque = df[(df.ID_ENTIDAD==9) & (df.ID_MUNICIPIO!=999)]
df_cdmx_parque = df_cdmx_parque.reset_index(drop=True)
Se realizó la suma de las columnas, con el objetivo de obtener todos los Autos, motos y camiones de pasajeros registrados en el parque vehicular.
df_cdmx_parque["AUTOSParque"] = df_cdmx_parque[['AUTO_OFICIAL','AUTO_PUBLICO','AUTO_PARTICULAR']].sum(axis=1)
df_cdmx_parque["MOTOSParque"] = df_cdmx_parque[['MOTO_OFICIAL','MOTO_DE_ALQUILER','MOTO_PARTICULAR']].sum(axis=1)
df_cdmx_parque["CAMParque"] = df_cdmx_parque[['CAM_PAS_OFICIAL','CAM_PAS_PUBLICO','CAM_PAS_PARTICULAR']].sum(axis=1)
df_cdmx_parque.head()
Se genero un nuevo dataframe para la posterior comparacion con los accidentes reportandos correspondientes a las categorias de autos, motos y camiones de pasajeros.
df_parque_vehicular = df_cdmx_parque[['ANIO','ID_MUNICIPIO','AUTOSParque','MOTOSParque','CAMParque']]
df_parque_vehicular.head(20)
Se realizó una agrupacion por años y municipios para posteriomente realizar una suma sobre las columnas del dataframe de accidentes para la posterior comparación.
df_accidentes = df_cdmx.groupby(["ANIO","ID_MUNICIPIO"])[['AUTOMOVIL', 'CAMPASAJ', 'MICROBUS', 'PASCAMION', 'OMNIBUS', 'TRANVIA',
'CAMIONETA','MOTOCICLET']].sum().reset_index()
Fue necesario sumar las columnas 'CAMPASAJ', 'MICROBUS', 'PASCAMION', 'OMNIBUS', 'TRANVIA', 'CAMIONETA' para la correcta comparacion del parque vehicular de los camiones de pasajeros.
df_accidentes["CAMIONPAS"] = df_accidentes[['CAMPASAJ', 'MICROBUS', 'PASCAMION', 'OMNIBUS', 'TRANVIA',
'CAMIONETA']].sum(axis = 1)
df_parque_accidentes = df_accidentes[['ANIO','ID_MUNICIPIO','AUTOMOVIL','MOTOCICLET','CAMIONPAS']]
df_parque_accidentes.head()
df_parque_accidentes = pd.merge(df_parque_accidentes,df_parque_vehicular)
df_parque_accidentes = df_parque_accidentes.rename(columns={'AUTOMOVIL': 'AUTOSAcci','MOTOCICLET':'MOTOSAcci','CAMIONPAS':'CAMAcci'})
df_parque_accidentes = pd.merge(df_parque_accidentes,delegacion_nombre[['ID_MUNICIPIO', 'NOM_MUNICIPIO']],how='left')
df_parque_accidentes.head()
Ya con el dataframe presente se comenzaron los análisis correspondientes.
Se empezo con el cálculo de estimados de locación y estadisticos de orden, obteniendo en la siguiente tabla el promedio, mediana, desviación estándar, min, max y cuartiles de los datos, con los datos sepaarados por alcaldia.
df_parque_accidentes[['AUTOSAcci','MOTOSAcci','CAMAcci','AUTOSParque','MOTOSParque','CAMParque']].describe()
Viendo estos números podemos inferir varias cosas:
1 En la columna de Accidentes de autos una gran parte de nuestros datos estan en valores menores a 1760.
2 En la columna de Accidentes de motos una gran parte de nuestros datos estan en valores menores a 85.
3 En la columna de Accidentes de camiones una gran parte de nuestros datos estan en valores menores a 436.
4 Se puede observar que los accidentes de todas las columnas la mayoria de los datos se concentran abajo de nuesto 75%.
5 En la parte del parque vehicular exiten valores muy grandes al llegar a nuestro maximo, ya que hay alcandias con una gran cantidad de vehiculos
En esta parte se agruparon los datos por año y de obtuvieron los mismos calculos del análisis anterior.
df_parque_accidentes.groupby('ANIO')[['AUTOSAcci','MOTOSAcci','CAMAcci','AUTOSParque','MOTOSParque','CAMParque']].sum().describe()
Tras estos cálculos se observa que la distribucion de los accidentes y parque vehicular por años esta mejor distribuido.
ax = sns.scatterplot(df_parque_accidentes['AUTOSParque'],df_parque_accidentes['AUTOSAcci']);
ax.set(xlabel = 'Parque vehicular Automoviles', ylabel = 'Número de accidentes')
ax.set_title('Gráfica análisis de Automovil Parque vs Accidentes')
ax = sns.scatterplot(df_parque_accidentes['MOTOSParque'],df_parque_accidentes['MOTOSAcci']);
ax.set(xlabel = 'Parque vehicular Motocicletas', ylabel = 'Número de accidentes')
ax.set_title('Gráfica análisis de Motocicleta Parque vs Accidentes')
ax = sns.scatterplot(df_parque_accidentes['CAMParque'],df_parque_accidentes['CAMAcci']);
ax.set(xlabel = 'Parque vehicular Camiones pasajeros', ylabel = 'Número de accidentes')
ax.set_title('Gráfica análisis de Camiones pasajeros Parque vs Accidentes')
Se observo que tanto en el scatter de motos y camiones exiten valores atípicos es por ello que se aplicaron filtros para posteiormente observar las respectivas correlaciones.
motos_filtro = df_parque_accidentes[df_parque_accidentes['MOTOSParque']>0]
ax = sns.scatterplot(motos_filtro['MOTOSParque'],df_parque_accidentes['MOTOSAcci']);
ax.set(xlabel = 'Parque vehicular Motocicletas', ylabel = 'Número de accidentes')
ax.set_title('Gráfica análisis de Motocicleta Parque vs Accidentes')
camiones_filtro = df_parque_accidentes[df_parque_accidentes['CAMParque']<6000]
ax = sns.scatterplot(camiones_filtro['CAMParque'],df_parque_accidentes['CAMAcci']);
ax.set(xlabel = 'Parque vehicular Camiones pasajeros', ylabel = 'Número de accidentes')
ax.set_title('Gráfica análisis de Camiones pasajeros Parque vs Accidentes')
Ya con las gráficas se obtuvieron las correlaciones entre las variables de parque vehicular y números de accidentes, obteniendo los siguientes resultados.
df_parque_accidentes['AUTOSParque'].corr(df_parque_accidentes['AUTOSAcci'])
motos_filtro['MOTOSParque'].corr(df_parque_accidentes['MOTOSAcci'])
camiones_filtro['CAMParque'].corr(df_parque_accidentes['CAMAcci'])
Se observa que entre el parque vehicular y el numero de accidentes exite una gran relación a partir del coeffiente obtenido, es por ello que se decide observar la relación entre todas las variables.
#Se eliminan las columnas innecesarias en nuestro análisis de correlaciones
filtro_motos = df_parque_accidentes['MOTOSParque']>0
df_total_filtro = df_parque_accidentes[filtro_motos]
df_total_filtro = df_total_filtro.reset_index(drop=True)
df_total_filtro_corr = df_total_filtro.drop(columns = ['ANIO','ID_MUNICIPIO','NOM_MUNICIPIO'])
df_total_filtro_corr.head()
df_total_filtro_corr.corr()
Con el fin de una mejor observación de la correlaciones entre las variables se hace uso de un heatmap.
plt.figure(figsize = (8,8))
ax = sns.heatmap(df_total_filtro_corr.corr(),annot=True, cmap = "YlGnBu",linewidths = 0.5)
Se observan correlaciones insospechadas, tales como la presente entre el número de accidentes en autos y en camiones de pasajeros, número de accidentes en autos y en motocicletas, etc., y se hizo uso de un pairplot para poder observar de mejor manera tales relaciones.
sns.pairplot(df_total_filtro_corr)
Debido a que el objetivo es demostrar la relación entre el parque y el numero de accidentes se deja para trabajos futuros el análisis de las relaciones entre accidentes.
from sklearn.linear_model import LinearRegression
lr = LinearRegression()
lr.fit(df_total_filtro['AUTOSParque'].to_frame(), df_total_filtro['AUTOSAcci'])
y_predict = lr.predict(df_total_filtro['AUTOSParque'].to_frame())
ax = sns.scatterplot(df_total_filtro['AUTOSParque'],y_predict);
sns.scatterplot(df_total_filtro['AUTOSParque'],df_total_filtro['AUTOSAcci'], ax=ax);
ax.set(xlabel = 'Parque vehicular Automoviles', ylabel = 'Número de accidentes');
ax.set_title('Gráfica regresión lineal de Automovil Parque vs Accidentes');
Podemos extraer la ecuación de la línea de nuestro modelo entrenado de esta manera. Y por lo tanto, nuestra función predictiva es la siguiente:
intercept = lr.intercept_
coefficient = lr.coef_[0]
print(f'y = {coefficient} * x + {intercept}')
efici = lr.score(df_total_filtro['AUTOSParque'].to_frame(), df_total_filtro['AUTOSAcci'])
print(f' Efectividad = {efici}')
lr.fit(df_total_filtro['MOTOSParque'].to_frame(), df_total_filtro['MOTOSAcci'])
y_predict = lr.predict(df_total_filtro['MOTOSParque'].to_frame())
ax = sns.scatterplot(df_total_filtro['MOTOSParque'],y_predict);
sns.scatterplot(df_total_filtro['MOTOSParque'],df_total_filtro['MOTOSAcci'],ax = ax);
ax.set(xlabel = 'Parque vehicular Motocicletas', ylabel = 'Número de accidentes')
ax.set_title('Gráfica regresión lineal de Motocicleta Parque vs Accidentes')
Podemos extraer la ecuación de la línea de nuestro modelo entrenado de esta manera. Y por lo tanto, nuestra función predictiva es la siguiente:
intercept = lr.intercept_
coefficient = lr.coef_[0]
print(f'y = {coefficient} * x + {intercept}')
efici = lr.score(df_total_filtro['MOTOSParque'].to_frame(), df_total_filtro['MOTOSAcci'])
print(f' Efectividad = {efici}')
lr.fit(df_total_filtro['CAMParque'].to_frame(), df_total_filtro['CAMAcci'])
y_predict = lr.predict(df_total_filtro['CAMParque'].to_frame())
ax = sns.scatterplot(df_total_filtro['CAMParque'],y_predict);
sns.scatterplot(df_total_filtro['CAMParque'],df_total_filtro['CAMAcci'],ax = ax);
ax.set(xlabel = 'Parque vehicular Camiones', ylabel = 'Número de accidentes')
ax.set_title('Gráfica regresión lineal de Camiones de pasajeros Parque vs Accidentes')
Podemos extraer la ecuación de la línea de nuestro modelo entrenado de esta manera. Y por lo tanto, nuestra función predictiva es la siguiente:
intercept = lr.intercept_
coefficient = lr.coef_[0]
print(f'y = {coefficient} * x + {intercept}')
efici = lr.score(df_total_filtro['CAMParque'].to_frame(), df_total_filtro['CAMAcci'])
print(f' Efectividad = {efici}')
Se observa que de entrre todas las regresiones la mas efectiva es la eslablecida a partir del análisis de los automoviles, además, que según las regresiones todos gráficas son con crecimiento, hay que tomar en cuenta qu este análsis se realizo tomando en cuenta todos los datos registrados por cada alcaldia.
El siguiente análisis se realizó tomando en cuenta los datos obtenidos unicamente por año, para observar la relación de una manera general.
df_total_anio = df_total_filtro.groupby('ANIO')[['AUTOSAcci','MOTOSAcci','CAMAcci','AUTOSParque','MOTOSParque','CAMParque']].sum()
df_total_anio.head()
lr.fit(df_total_anio['AUTOSParque'].to_frame(), df_total_anio['AUTOSAcci'])
y_predict = lr.predict(df_total_anio['AUTOSParque'].to_frame())
ax = sns.scatterplot(df_total_anio['AUTOSParque'],df_total_anio['AUTOSAcci']);
sns.scatterplot(df_total_anio['AUTOSParque'],y_predict,s=40,ax = ax);
ax.set(xlabel = 'Parque vehicular Automoviles', ylabel = 'Número de accidentes');
ax.set_title('Gráfica regresión lineal de Automovil Parque vs Accidentes anual');
Se observa una tendencia a la baja en relacion al parque vehicular analizando los datos de manera anual.
lr.fit(df_total_anio['MOTOSParque'].to_frame(), df_total_anio['MOTOSAcci'])
y_predict = lr.predict(df_total_anio['MOTOSParque'].to_frame())
ax = sns.scatterplot(df_total_anio['MOTOSParque'],df_total_anio['MOTOSAcci']);
sns.scatterplot(df_total_anio['MOTOSParque'],y_predict,s=40, ax= ax);
ax.set(xlabel = 'Parque vehicular Motocicletas', ylabel = 'Número de accidentes')
ax.set_title('Gráfica regresión lineal de Motocicleta Parque vs Accidentes anual')
Se observa una tendencia a la alza en relacion al parque vehicular analizando los datos de manera anual.
lr.fit(df_total_anio['CAMParque'].to_frame(), df_total_anio['CAMAcci'])
y_predict = lr.predict(df_total_anio['CAMParque'].to_frame())
ax = sns.scatterplot(df_total_anio['CAMParque'],df_total_anio['CAMAcci']);
sns.scatterplot(df_total_anio['CAMParque'],y_predict,s=40,ax = ax);
ax.set(xlabel = 'Parque vehicular Camiones', ylabel = 'Número de accidentes')
ax.set_title('Gráfica regresión lineal de Camiones de pasajeros Parque vs Accidentes anual')
Se observa una tendencia a la baja en relacion al parque vehicular analizando los datos de manera anual.
Por lo observado anteriormente se decidio realizar visulaciones de las tres alcaldias con mayor número de accidentes reportados.
Debido a que el objetivo es demostrar la relación entre el parque y el numero de accidentes se deja para trabajos futuros el análisis de las relaciones entre accidentes.
alcaldia1 = df_total_filtro['NOM_MUNICIPIO']=='Gustavo A. Madero'
alcaldia2 = df_total_filtro['NOM_MUNICIPIO']=='Iztapalapa'
alcaldia3 = df_total_filtro['NOM_MUNICIPIO']=='Cuauhtémoc'
df_total_alcaldia = df_total_filtro[(alcaldia1) | (alcaldia2) | (alcaldia3)].reset_index()
ax = sns.scatterplot(df_total_alcaldia['AUTOSParque'], df_total_alcaldia['AUTOSAcci'], hue=df_total_alcaldia['NOM_MUNICIPIO']);
ax.set(xlabel = 'Parque vehicular Automoviles', ylabel = 'Número de accidentes');
ax.set_title('Gráfica Automovil Parque vs Accidentes');
ax = sns.scatterplot(df_total_alcaldia['MOTOSParque'], df_total_alcaldia['MOTOSAcci'], hue=df_total_alcaldia['NOM_MUNICIPIO']);
ax.set(xlabel = 'Parque vehicular Motocicletas', ylabel = 'Número de accidentes');
ax.set_title('Gráfica de Motocicleta Parque vs Accidentes');
ax = sns.scatterplot(df_total_alcaldia['CAMParque'], df_total_alcaldia['CAMAcci'], hue=df_total_alcaldia['NOM_MUNICIPIO']);
ax.set(xlabel = 'Parque vehicular Camiones', ylabel = 'Número de accidentes');
ax.set_title('Gráfica de Camiones de pasajeros Parque vs Accidentes');
Se observa que en la parte de automovil las tres alcaldias presentan comportamientos similares tendiendo a la baja, mientras que en la parte de motocicletas la tendencia en un poco parecida entre ellas y por último en las observaciones de los camiones de pasajeros sus comportamietnos son parecidos a pesar de la diferencia entre el parque vehicular.
Se agregan las gráficas de los comportamientos de cada tipo de vehiculo por alcaldia para que el lector pueda observar con se comportan individualmente.
fig, axes = plt.subplots(4, 4, figsize=(25, 15));
indice = 0;
for i in range(0,4):
for j in range(0,4):
alcaldia = df_total_filtro['NOM_MUNICIPIO']==df_total_filtro['NOM_MUNICIPIO'].unique()[indice]
sns.scatterplot(df_total_filtro[alcaldia]['AUTOSParque'],df_total_filtro[alcaldia]['AUTOSAcci'],ax=axes[i, j])
axes[i, j].set(xlabel='', ylabel='', title=df_total_filtro['NOM_MUNICIPIO'].unique()[indice]);
indice = indice+1
fig.suptitle('Graficas Automovil', fontsize=25);
fig, axes = plt.subplots(4, 4, figsize=(25, 15));
indice = 0;
for i in range(0,4):
for j in range(0,4):
alcaldia = df_total_filtro['NOM_MUNICIPIO']==df_total_filtro['NOM_MUNICIPIO'].unique()[indice]
sns.scatterplot(df_total_filtro[alcaldia]['MOTOSParque'],df_total_filtro[alcaldia]['MOTOSAcci'],ax=axes[i, j]);
axes[i, j].set(xlabel='', ylabel='', title=df_total_filtro['NOM_MUNICIPIO'].unique()[indice]);
indice = indice+1
fig.suptitle('Graficas Motocicletas', fontsize=25);
fig, axes = plt.subplots(4, 4, figsize=(25, 15));
indice = 0;
for i in range(0,4):
for j in range(0,4):
alcaldia = df_total_filtro['NOM_MUNICIPIO']==df_total_filtro['NOM_MUNICIPIO'].unique()[indice]
sns.scatterplot(df_total_filtro[alcaldia]['CAMParque'],df_total_filtro[alcaldia]['CAMAcci'],ax=axes[i, j]);
axes[i, j].set(xlabel='', ylabel='', title=df_total_filtro['NOM_MUNICIPIO'].unique()[indice]);
indice = indice+1
fig.suptitle('Graficas Camiones de pasajeros', fontsize=25);
Para esta última parte y como se resalto en múltiples comentarios los análisis inferenciales deberán ir enfocados en evaluar 3 aspectos específicos.
El primer rubro será el análisis del impacto de los distintos factores o variables sobre los choques, para esto será necesario llevar a cabo análisis de regresión multivariable en los cuáles se puedan evaluar los coeficientes, así como las interacciones de estos factores. Dentro de las herramientas estadísticas que podremos utilizar para este análisis destacan ANOVA, K - Means, Random Forests, entre otros.
El segundo rubro será ahondar en el análisis comparativo con otros archivos de datos como lo es el parque vehicular, la cantidad de hombres y mujeres registrados, persona registradas por delegación, personas registradas por grupos de edad, entre otros. Los cuáles serán eficacez para eliminar ciertos tipos de sesgo que por ahora se tienen y que impiden un análisis más profundo de los factores determinates para la generación de un choque.
Como último rubro, será importante seguir enfocando hacía la predicción de choques bajo diferentes tipos de escenarios, para eto será útil efectuar métodos como las series de tiempo. asimismo y complementadno lo anterior estas predicciones podrán ayudar a también generar posibles análisis económetricos de la situación.
Finalmente, englobando todo este futuro desarrollo, estipulado sobre las base del análisis exploratorio que se genero, la propuesta final será diseñar estrategias que permitan mitigar los factores (causas raíz) que propician mayor número de choques. A manera de ejemplo, durante este primer análisis encontramos que a medida que incrementaba la hora del día, aumentaba el número de choques, y que el Viernes era el día con mayor número de choques registrados. Aunado a esto encontramos, qué alcaldías, erán las que tenían un reporte de mayor número de accidentes, por lo que, considerando estos tres factores, podriamos diseñar una medida de contigencía de alchólimetros mejor ubicados y distribuidos sobre las delegaciones con más concentración de choques y en horarios que quizas antes no se tenían contemplados, con esto probablemente reducir choques y consecuentemente reducir gasto público e indivdual. Esto claramente fue un somero esbozo del tipo de soluciones que podrían proponerse, siendo conscientes de que ese tipo de decisiones deberá contemplar también otros factores como movilidad en la zona, presupuesto de tránsito y vialidad, entre otros que los expertos en el tema considerarán, sin embargo, este ejercicio tenía como mera finalidad evidenciar el potencial de los resultados que pudiera llegar a tener este análisis y la utilidad sobre la problemática planteada al inicio del proyecto.
Modelo ARIMA
ARIMA es un modelo autorregresivo integrado de promedio móvil (acrónimo del inglés autoregressive integrated moving average) este modelo utiliza variaciones y regresiones de datos estadísticos con el fin de encontrar patrones para lograr efectuar una predicción futura del comportamiento de la serie.
En este caso específico se realizo un modelo ARIMA con el fin de evaluar el comportamiento del número de accidentes a lo largo de los años venideros, para ello se realizaron varias modificaciones a nuestro Data Frame principal "df_cdmx", a continuación se presenta la descripción detallada del desarrollo de este modelo.
Lo primero fue agrupar el número de accidentes por año y mes, para obtener una suficiente cantidad de datos que nos permitiera llevar a cabo el análisis, esto debido a que el análisis por años tenía unicamente 10 registros.
filtro_arima=df_cdmx.groupby("ANIO")["MES"].value_counts()
filtro_arima1=filtro_arima.reset_index(name="Numero_Acc")
filtro_arima2 = pd.DataFrame(filtro_arima1)
filtro_arima2
Lo siguiente fue establecer como índice a las columnas de Año y Mes, asimismo ordenamos los meses (1-12) por año para proceder con la graficación del DF con el que estariamos trabajando.
filtro_arima3 = filtro_arima2.set_index(["ANIO","MES"])
filtro_arima3_1 = filtro_arima2.set_index(["ANIO"])
filtro_arima3.head(12)
filtro_arima4 = filtro_arima3.sort_values(by=['ANIO','MES'])
filtro_arima4_1 = filtro_arima3_1.sort_values(by=['ANIO','MES'])
filtro_arima5 = pd.DataFrame(filtro_arima4)
filtro_arima5
fig, ax = plt.subplots()
filtro_arima5.plot(ax=ax)
plt.show;
Lo siguiente fue llevar a cabo una análisis para identificar si necesitamos transformar nuestra Xt esto con el fin de eliminar la no estacionalidad en la media p en la heteroscedasticidad. Para ello utilizamos el Dicky-Füller, estableciendo una hipótesis que de no ser rechazada permitía continuar con el proceso sin necesidad de realizar el ajuste anteriormente mencionado. (P-value > .05 por lo tanto continuamos)
from statsmodels.tsa.stattools import adfuller
result = adfuller(filtro_arima5['Numero_Acc'])
print('ADF Statistic: {}'.format(result[0]))
print('p-value: {}'.format(result[1]))
print('Critical Values:')
for key, value in result[4].items():
print('\t{}: {}'.format(key, value))
Definimos nuestro Train y Test set
train = filtro_arima5.iloc[0:100]
test = filtro_arima5.iloc[101:]
Utilizamos ARIMA del módulo statsmodels para evaluar el modelo de ajuste, el método, el análisis logarítmico y los componentes con sus respectivos coeficientes de autoregresión y promedios móviles
from statsmodels.tsa.arima_model import ARIMA
modelo=ARIMA(train['Numero_Acc'],order=(1,0,5))
modelo=modelo.fit()
modelo.summary()
inicio =len(train)
final = len(train) - len(test)
pred = modelo.predict(inicio = inicio , final = final,typ='levels').rename('ARIMA Predictions')
pred.plot(legend=True);
test['Numero_Acc'].plot(legend=True)
Pese a que en el análisis decriptivo del modelo de ARIMA en donde se habían identificado la raíces MA.1 - MA.5 (Moving Averages), con coeficientes ue se muestran en el indicador como significativos el AR.1 pudiera ser indientificado como no significativo con un valor de z de 35, en este caso resultaría efectivo analizar de nuevo nuestros sets de entrenamiento y prueba y el modelo ARMA (1,5).